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ABSTRACT 

Massive stars provide feedback that shapes the interstellar medium of galaxies at 
all redshifts and their resulting stellar populations. Here we present three adaptive 
mesh refinement radiation hydrodynamics simulations that illustrate the impact of 
momentum transfer from ionising radiation to the absorbing gas on star formation 
in high-redshift dwarf galaxies. Momentum transfer is calculated by solving the ra- 
diative transfer equation with a ray tracing algorithm that is adaptive in spatial and 
angular coordinates. We find that momentum input partially affects star formation 
by increasing the turbulent support to a three-dimensional rms velocity equal to the 
circular velocity of early haloes. Compared to a calculation that neglects radiation 
pressure, the star formation rate is decreased by a factor of five to 1.8 x 10~ 2 M Q 
yr _1 in a dwarf galaxy with a dark matter and stellar mass of 2.0 x 10 8 Af© and 
4.5 x 10 5 M Q , respectively, when radiation pressure is included. Its mean metallicity 
of 10 -2 ' 1 Zq is consistent with the observed dwarf galaxy luminosity-metallicity re- 
lation. However, what one may naively expect from the calculation without radiation 
pressure, the central region of the galaxy overcools and produces a compact, metal- 
rich stellar population with an average metallicity of 0.3 Z Q , indicative of an incorrect 
physical recipe. In addition to photo-heating in H II regions, radiation pressure further 
drives dense gas from star forming regions, so supernovae feedback occurs in a warmer 
and more diffuse medium, launching metal-rich outflows. Capturing this aspect and a 
temporal separation between the start of radiative and supernova feedback are numer- 
ically important in the modeling of galaxies to avoid the "overcooling problem" . We 
estimate that dust in early low-mass galaxies is unlikely to aid in momentum transfer 
from radiation to the gas. 

Key words: cosmology — methods: numerical — hydrodynamics — radiative trans- 
fer — star formation 



1 INTRODUCTION 

Stellar radiation from massive stars and their supernova 
(SN) explosions significantly alter the surrounding interstel- 
lar medium (ISM) and, thus, sub sequent star formation i n 
the host galaxy (for a review, see iMcKee fc OstrikenEfflOTh . 
On a grander scale, their input can drive large-scale outflows 
that increas e the entropy and metallicity of the intergalac- 
tic m edium (|Davies et al.ll998l ; lAguirre et al j2001al ;lB enson 
I2010I . for a review). To date, the primary feedback mecha- 
nism in galaxy formation simulations originate from SNe, 
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which is implemented by locally injecting ther mal and ki- 
netic energy and metals |Cen fc Ostrikerl 19921 ). T here are 
vario us relevant physics, such as photoheat ing (e.g. iGnedinl 
l200Cf ). momentu m input from radiation fe.g. lHaehneltll995l ). 
cosmi c rays (e.g. lJubelgas et al.ll200sl) and magnetic fields 
(e.g. lWangfcA bcl 2009; Kotarba et"al]|201ll 'l. that can play 
an important role in star and galaxy formation. Numerical 
simulations are a useful tool in exploring the impacts of these 
physical processes in galaxy formation. But only recently, 
feedback mechanisms other than SNe have begun to appear 
in cosmological simulations, sometimes in a phenomenolog- 
ical manner, because of resolution and algorithmic limita- 
tions. 
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Simulating the first low-mass galaxies with the cor- 
rect physical model is crucial because they are the build- 
ing blocks for all subsequent galaxies, and it is clear that 
stellar radiation from them were pr e dominately responsible 
for r e ionisation (e.g. ISnapird 1 19861 : iFan. Carilli fe Keating! 



120061 : n ouwens et al.l l201ll ) . The most luminous of these 
have been detected at 6 < z < 10 in the 



reiomscrs 



campaigns (e.g. 


Bouwens et all 20081: McLure et al.l l20ld: 


iFinkelstein et al. 


20ld: Ouchi et all 201Ch. Even lower mass 



and higher redshift galaxies will be detected with the At- 
acama Large Millimeter Array (ALMA) and the James 
Webb Space Telescope (JWST), providing additional con- 
straints on early galaxy formation, reionisation, and the 
early metal enrichment of the intergalactic medium. Closer 
to home, nearby dwarf galaxies provide clues to their for- 
mation, , where_jdl_of_thj^^ or before reioni- 
sation (|Grebel. Gallagher fc Harbeckll2003l). in the ir [a/Fe] 
versus [Fe/Hl evolution jTolstov. Hill fc Tosil 120091 ). metal- 
licity distributio n functions (|Kirbv et al.l 1201 lh . metallicitv 
radia l gradients ( Tolstoy et al.l 20041 ) . star formation histo- 



ries llMonelli et all |2010|) . and individual ab undance pat- 



terns |SebdeTai]l2oTa; ICaffau et al1l201ll ). These dat 



from the high-redshift and the local universe provide im- 
portant constraints on numerical models of star formation 
and feedback. 

The first generation of gal axies with 10 8 < M^aio/ 'Mq < 
10 10 (|Bromm fc Yoshidall201ll) are an excellent cosmological 
testbed to numerically investigate the impact of each phys- 
ical process. Compared to simulating larger galaxies, they 
present a relatively clean scenario with few instances of prior 
star formation that sets the "initial conditions" for galaxy 
formation. They are also small enough so it is computation- 
ally feasible to achieve high enough resolution to well resolve 
star forming regions. Their progenitors may have hosted sev- 
eral metal-free and massive (Pop III) stars, enriching the 
IGM and pre-heating the gas to thousands of degrees, which 
then re-accret es into the dark matter (DM) halo and forms a 
dwarf galaxy l|Wise fc Abelll200Sl : iGreif et al.ll2010l ). In this 
paper, we investigate the impact of radiative cooling from 
metals and a Hb-dissociating radiation background, but the 
main focus of this work is the effect of momentum trans- 
fer from ionising radiation, i.e. radiation pressure, on star 
formation in such galaxies. 

Regulating star formation by radiation pressure is not 
a new idea, esp ecially in the Mi l ky Way and present-day 
star formation dShu et~aU Il99ll : iFerraral fl99a: iLi fc Shul 
1 19961 : iMatzner fc McKeeHl999l ). However, lHaehneltl (|l995l ) 
was the first to apply this idea to galaxy formation. He 
found that radiation pressure could be the dominant feed- 
back process within ~110 pc in DM haloes with masses 
< 10 9 — lO lo A/0. He suggested that it could be more im- 
portant than SNe heating in driving outflows from such 
early dwarf galaxies. In turn, this would lower their gas 
fractions and alter their predicted mass-to-light ratios. In 
larger galaxies, it was unclear whether momentum deposi- 
tion could drive outflows through ionising radiation alone. 
However, UV radiation can be absorbed by dust and be 
re-emitted many times, increasing the momentum trans- 
fer to the absorber. Without making the assumption that 
gas dynamics are coupled to dust grains, radiation pres- 
sure from massive stars were found to drive dust outflows 



in cosmological simulations, enriching the IGM to the C 
and Si abundances fo und i n z = 3 quasar absorption lines 
llAguirre et al. _ 200 1 a , b | c L_ Murray. Quataert fc Thompson! 



i|2005l ) and iThompson, Quataert fc Murray! (|2005h investi- 
gated the effects of radiation pressure in driving galactic out- 
flows from starburst galaxies and regulating star formation 
by reducing the gas supply for further star formation and 
BH growth. They found that it could explain the Mbh — cr 
relation for early-type galaxies, where a is the stellar velocit> 



lty 

dispersion dFerrarese fc MerrittlbOOol : iGebhardt et al.ll200oi : 



iTremaine et all 12002 ) . They also found that momentum- 
driven winds could be launched if the opacity k > 6 cm 2 g _1 
in a starburst galaxy, assuming that momentum transfer is 
increased by a factor of Tin — 1, where rm is the optical 
depth in the infrared (IR). 

These works led to various implementations of mo- 
mentum transfer in numerical studies of galaxy for- 
mation, usually in the form of "kicking" particles 
with some wind velocity f ro m star-forming region s (e.g . 
Springel fc Hernquistl 120031: lOppenheimer fc Pavel 2008 



Sales et al.l 120101 : iHopkins. Quataert fc Mur ray 2011). The 



early implementations decoupled these wind particles from 
the hydrodynamics to create outflo ws; otherwise, they would 
be con fined to the galaxy. Recently. iKrumholz fc Thompson] 
|2012l) investigated the parsec-scale effects of momentum 
transfer in radiation hydrodynamics simulations and com- 
pared them to the aforementioned subgrid models. They 
found that most of the radiation is not trapped because the 
optically-thick shell breaks up from a radiation Rayleigh- 
Taylor instability, i.e. momentum transfer does not scale as 
rm- This may lead to an overestimate of momentum trans- 
fer in the subgrid models when the radiation source shines 
above 10 per cent of its Eddington limit. 

Here we present results that demonstrate the effects 
of momentum transfer from ionising radiation on dwarf 
galaxy formation, using cosmological radiation hydrody- 
namics simulations. The amount of momentum deposition 
is directly calc ulated from our ra diation transport solver 
enzo+moray (|Wise fc Abelll201ll ). In the next section, we 
outline some analytics that govern the dynamics of outflows 
and turbulence created by radiation pressure. Then in Sec. [31 
we describe the numerical simulation setup and algorithms. 
In Sec. [4] we present the effects of radiation pressure on star 
formation and metal enrichment in dwarf galaxies, prevent- 
ing the overcooling problem seen often in galaxy formation 
simulations. We discuss possible constraints on high-redshift 
dwarf galaxies and then focus on differences in feedback im- 
plementations in simulations in Sec. [5] We finish that section 
with an estimate of momentum transfer onto dust grains in 
dwarf galaxies. We summarise our findings in SecO 



2 MOMENTUM TRANSFER FROM IONISING 
RADIATION 

First we review some basic analytical arguments 
on momentum transfer from ionising r a diation . We 
similar approach a s Haehneltl |l995) and 



take a , , 

iMurrav. Quataert fc Thompson! 1 20051 ) and use an 
Eddington-like argument on an entire galaxy, equating 
the radiation pressure from ionising stellar radiation to 
the gravitational potential of the host dark matter halo. 
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Its application to local sta r formation and feedbac k has 
been extensively studied (|McKee fc Qstrikerl 120071 . and 
references therein). But can we apply the same principles 
to an entire dwarf galaxy? How will the energetics differ if 
the momentum transfer originates from ionising radiation 

instead of protostellar winds or Thomson scatte ring? 

In the traditional Eddington approach (Ed dingtonl 

1 1926h . Thomson scattering is the source of radiation pres- 
sure. The Eddington luminosity limit is 

4irGMm p c 

^edd — , (1) 

where M is the mass of object, m p is the proton mass, 
and <tt = 6.65 x 10~ 25 cm 2 is the Thomson cross-section. 
However, the absorption cross-section for hydrogen <thi = 
6.35 x 10 -18 cm 2 at 1 Ryd, which is ~ 10 7 times larger 
than the Thomson cross-section. However, it should be noted 
that a spectrum-weighted cross-section is somewhat lower 
because cthi oc t>~ 3 . Thus momentum transfer during the 
absorption of ionising photons might significantly affect the 
large-scale dynamics in the galaxy. 

For illustrative purposes, we focus on a single expanding 
H II region with hydrogen only that is embedded in a DM 
halo. A shock front forms at the edge of the H H region 
as w hen the ionisation front becomes D-type (|Osterbrockl 
1989), i.e., the expansion speed of the H II region slows below 
the sound speed of the ionised medium. If the ionisation 
front coincides with the shock front, the optical depth from 
the radiation source will exceed unity within the shock. We 
can simplify the calculation by assuming that the H II region 
is optically thin and all of the ionising radiation momentum 
is tr ansferred to the optically-thick shock. 

iMurrav. Quataert fc Thompson! (|2005h considered the 
spherically symmetric and self-gravitating case where the 
stellar system is contained within a halo with a gas frac- 
tion /gas- Their main focus was on nuclear starbursts and 
the resul ting outflows and explaining the M — a rela- 
tionship jFerrarese fc Merrittl l200fj| ; iGebhardt et all l200d ; 



i Tremaine ' et al l I 2 002^ ~rle7 refected thermal pre™ 



forces from the H 11 region because it is negligible for a stellar 
system that is centered in the gravitational potential of halos 
that host L > L* galaxies. For a central ionising luminos- 
ity L where each photon is absorbed once, the momentum 
equation is 



dP 
~dt 



GM(r)M gas (r) 



+ 



(2) 



In the case where photons are absorbed and re-emitted many 
times in a dusty medium, the second ter m can become im- 
portant. IMurrav. Quataert fc Thompson! showed that for a 
luminosity L > 4/ gas a 4 c/G, the gas moves outward driven 
by pure radiation pressure, where a is the halo velocity dis- 
persion and all radiation is assumed to be absorbed. 

For a dwarf galaxy with a = 20 km s _1 and the cos- 
mic baryon fraction, this critical luminosity is on the or- 
der of 10 42 erg s" 1 (1O 84 L ), equivalent to a mass-to-light 
ratio on the order of unity for a Salpeter IMF. The grav- 
itational potential in dwarf galaxies with V c < 30 km s" 1 
cannot prevent the expansion of the dense shell surround- 
ing H 11 regions into the IGM. Furthermore, star formation 
is not necessarily concentrated in the centre of these small 
irregular galaxies jTolstov. Hill fc Tosil200gh . If th e H 11 re- 
gion is pressure-supported (e.g.. see lFerrars 1 1993! ). then we 
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Figure 1. Radius (top) and velocity (bottom) of an expanding 
optically thick shell driven by momentum transfer from ionising 
radiation for a central luminosity of L = 10 5 , 10 6 , 10 7 Lq with an 
average ionising photon energy of 20 eV. 



can neglect the gravitational force term in Eq. [2] Taking 
L to be time-independent and instantaneously non-zero at 
t > 0, the momentum equation for the optically-thick shell is 
simply P(t) = Lt/c. We consider the spherically-symmetric 
case where the shell sweeps all of the ISM with a density 
p = p(r) in its path, and it moves at a velocity 



«(*) = 



Lt 



cM s] 



Lt 

c 



4tt 



r' 2 p(r')dr' 



(3) 



Solving for radius and velocity with a homogeneous medium 
p(r) = po and isothermal density profile p(r) = po(r /ro)" 2 , 
we obtain 

,1/2 



r(t) 



(rf + 2At^ 
(r 2 + At 2 /3r 2 ) J 



(P = Po) 
(p oc r~ 2 ) 



(4) 



«(*) = 



tA (rf + 2At 2 ) 
t(A/3r 2 ) (r 2 + ^ 2 /3r 2 ) 



-3/4 
-1/2 



= Po) 
(p oc r -2 ) 



(5) 



where A = 3L / 4n poc, and r; is the radius where the optically 
thick shell first forms. The initial Stromgren sphere radius is 
calculated by equating the ionising photon luminosity iV 7 = 
L/{E y ) to the number of recombinations, so that 

1/3 

(6) 



/ 3iV 7 



given an initial absorber number density no of the star- 
forming region. Here (-E 7 ) is the average ionising photon 
energy, and qb » 2.5 x 10 -13 cm 3 s _1 is the Case B recom- 
bination rate at T = 10 4 K. Molecular complexes are very 
clumpy and, in the Galactic plane, they have an average 
densi ty in the range of 10 2 — 10 4 cm" " 3 jBergin fc Tafallal 
2007). For a representative stellar cluster, 



n = 1.57 Lg /3 

where = 

n3 -3 



1/3 -1/3 



-^7,20 n 3 

L/10 6 L 



0> 



pc, 

-£7,20 



(7) 



£ 7 /20 eV, and n 3 



10 5 ,10 6 ,10 7 L Q , 



no/10 3 cm" 3 . We plot r(t) and v(t) for L 
po = 1 cm" 3 , and (£7) = 20 eV in Figure [T] for a uni- 
form and a isothermal density profile. We set rg = 15 pc 
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so that both the uniform and isothermal sphere with ra- 
dius 3ro contains 10 4 Mq. In a uniform medium, the shell 
rapidly accelerates to its maximum velocity ~ 100 km s _1 af- 
ter 10 4 yr. Its velocity decreases as the shell accumulates 
material from the ISM. It eventually reaches 1 kpc after 
1 Gyr for the 10 6 Lq case. In the isothermal sphere, the 
shell accelerates as it runs down the density slope, reaching 
1 kpc after 50 Myr. The shell reaches a constant velocity 
v 3 = (L/iiTpqcrn) 1 ^ 2 ~ 13L6 km s~* after about a Myr 
( Bally fc Ladal 1 19831 ; IShu et all Il99ll ; iMatzner fc McKeei 



3 COSMOLOGICAL SIMULATIONS 



L 6 0,5 of the radiation 



Il999h with"afraction / = 5.0 x 10 
energy being transferred to the kinetic energy of the shell. 

This ideal model will break down as the shell when 
gravitational forces from the DM halo becomes larger than 
L/c and the thermal pressure inside the H II region. How- 
ever if parts of the shell fragment because of ionisation 
front or thin-shell instabilities, the diffuse gas in the open- 
ings of the shell will experience a greater acceleration, i.e. 
cham pagne flows (e.g. lYorke. Tenorio-Tagle fc Bodenheimer] 
Il983h . These radiation pressure driven champagne flows are 
a possible mechanism for enhancing mass outflow velocities 
in dwarf galaxies in this model. 

The terminal velocity for the shell traveling down a r~ 2 
density profile is comparable to the circular velocity of the 
DM halo. Radiative pressure driven winds in dwarf galaxies 
may contribute more to the driving of turbulence instead of 
outflows when the shells are contained within the halo. Colli- 
sions between surviving swept-up shells can drive turbulence 
in the star-forming region. This has b een shown to be like! 



3 

in present-d a y star formation (e.g. iNorman fc Silk 1980; 
Franco! [l983l ; iKrumholz. Matzner fc McKeei 120061 ; iMatznerl 



20071 ). This added turbulent support may help in regulating 



star formation by preventing a catastrophic central collapse 
in a galaxy. 

We can understand the evolution of turbulence by in- 
vestigating the vorticity equation, 



~Dt 



-wV ■ v 



VPxVp 2 
h vv u, 



(8) 



where D/Dt is the fluid derivative, V is pressure, and v is 
the viscosity. In our analysis, we consider non-viscous fluids 
and set v = 0. The first term describes the stretching of 
vorticity; the second term comes from non-barotropic flows, 
i.e. V 7^ V{p), and the last term accounts for the dissipa- 
tion of turbulence through viscous forces. Non-barotropic 
flows occur at and near shock fronts and have been found to 
be the main process that sustains vorticity in several situ- 
ations. For example, virial shocks and internal weak shocks 
have been found to drive turbulence dur ing halo virialization 
jWise fc Abelll2007al ; iGreif et all 120081 '). Furthermore, non- 
thermal pressure, such as radiation pressure or cosmic rays, 
contributes to the generation of vorticity when its gradient 
(acceleration) is not aligned with the density gradient. To 
determine the importance of radiation pressure in outflows 
and sustaining turbulent motions, we use cosmological radi- 
ation hydrodynamics simulations that focus on dwarf galaxy 
formation, which we describe next. 



We ex pand upon the simulation presented in W iseetajj 
l|2012l . hereafter Paper I) that focused on the transition from 
Population III to II star formation and the effects of radia- 
tive and supernova feedback on early dwarf galaxy forma- 
tion. In this work, we compare two new simulations that add 
radiative cooling from metals, an Eb-dissociating radiation 
background, and momentum transfer from ionising radia- 
tion to the original physics set. These simulations were run 
with the adaptive mesh r efinement (AMR) c ode enzo v2.C0 
ijBrvan fc Normanlll997l : lO'Shea et aT1l2004l). It uses an TV- 
body adap t ive pa r ticle-mesh solver | Efstathiou et al.l 1 19851 ; 
ICouchmanl Il99ll ; iBrvan fc Normanl I1997m ~ to follow the 
DM dynamics. It solves the hydrodynamical equation us- 
ing the second-order accurate piecewise para bolic method 
jWoodward fc Colellal[l98i : lBrvan et alJll995j ). while a Rie- 
mann solver ensures accurate shock capturing with mini- 
mal vi scosity. We use the recently a dded HLLC Riemann 
solver (|Toro. Spruce fc Speareslll994l ') for additional stabil- 
ity in strong shocks and rarefaction waves. 

We use a simulation box of 1 Mpc on a side with a base 
resolution of 256 3 , resulting in a DM mass resolution of 1840 
M Q . This allows us to resolve 2 x 10 5 M S DM haloes t hat 
can cool by H2 in the gas phase dTegmark et al.lll997h by 
at least 100 particles. We refine the grid on baryon overden- 
sities of 3 x 2 _0,2i , where I is the AM R level, resulting in 
a sup er-Lagrangian behaviour (also see lO'Shea fc Normanl 
200S)- We also refine on a DM overdensity of three and 
always resolve the local Jeans length by at least 4 cells, 
avoiding artific i al fra gmentation during gaseous collapses 
jTruelove et al.lll997l ). If any of these criteria are met in 
a single cell, it is flagged for further spatial re finement. 

W e initialise the simulations with GRAFIC (|Bertschingerl 
l200lft at z = 130 and use the cosmological parame- 
ters from the 7-year WMAP ACDM+SZ+LENS best fit 
|Komatsu et al.l 120111 ): Q M = 0.266, n A = 0.734, ft 6 = 
0.0449, h = 0.71, a 8 = 0.81, and n = 0.963 with the vari- 
ables having their usual definitions. We use a maximum re- 
finement level of I = 12, resulting in a maximal comoving 
resolution of 1 pc. We do not smooth the DM mass field 
at any level, which is done to reduce artificial gas heating 
from DM particles in si mulations that have regions th at are 
baryon dominated (e.g. lAbel. Bryan fc Normanl [2002'). The 
adaptive particle-mesh solver has a force resolution of two 
cell widths of a given AMR grid. We compare the simula- 
tions at z = 8.1, which is when the metal cooling simulation 
reaches an ionised volume fraction of 96 percent. 



3.1 Star formation 

We distinguish Pop II and Pop III SF by the total metallic- 
ity of the densest cell in the molecular cloud. Pop II stars 
are formed if [Z/H] > —4, and Pop III stars are formed oth- 
erwise. We do not consider hypothetical Pop III. 2 stars or 
intermediate mass stars from CMB-limited cooling. Here we 
outline how we treat Pop II and III SF in our simulations; 
for more details and justifications, we refer the reader to 
IWise et al.l l|2012l ). 



1 enzo.googlecode.com, changeset 03a72f 4f 189a 
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3.1.1 Population III star formation 

We use the same Pop III SF model as IWise fc Abel (|2008l ) 
where each star particle represents a single star. In this 
model, a star particle forms when a cell has all of the fol- 
lowing criteria: 

(i) an overdensity of 5 x 10 5 (~ 10 3 cm" 3 at z = 10), 

(ii) a converging gas flow (V • v gas < 0), and 

(iii) a molecular hydrogen fraction /h2 > 5 x 10 -4 . 

These physical conditions are typical of collapsing metal- 
free molecular clouds ~ 10 Myr before the birth of a 
Pop III main-sequence s tar l)Abel. Brvan fc Normanl 120021 : 
lO'Shea fc Normanll2007l ). This is approximately the end of 
the initial free-fall phase, the so-called "loitering phase", 
where the minimum temperature a t this point characterise s 
the fragmentation mass scale (e.g. lOmukai fc Yoshii| [2003). 
This prescription is similar to the Cen fc Qstrikerl ( 19921 ) 
method but removes the Jeans unstable (M gas > Mj) and 
cooling timescale (i coo i < idyn) requirements. We do not 
consider the former criterion because it is not applicable to 
simulations that resolve the Jeans length at all times. The 
molecular hydrogen fraction requirement effectively con- 
strains star formation to cooling molecular clouds, where 
the H2 formation rate is significantly larger than the disso- 
ciation rate from a Lyman- Werner radiation field. 

If multiple cells meet the star particle formation criteria 
within 1 pc, we form one Pop III star particle with the centre 
of mass of these flagged cells to ensure that one massive 
star is created per metal-free molecular cloud. We randomly 
sample from an IMF with a functional form of 



/(log M)dM = M~ 



exp 



/ Mchar V 
KM) 



dM 



(9) 



to determine the stellar mass. Above Af c har, it behaves as 
a Salpeter IMF but is exponentially cutoff below that mass 
l|Chabrierll2003l ). After the star particle forms and its mass 
is determined, an equal amount of gas is removed from the 
computational grid in a sphere that contains twice the stel- 
lar mass and is centred on the star particle. At this stage 
of the collapse, a sphere enclosing 200Mq has a radius of 
1-2 pc and has approximately 200 cells. The star particle 
acquires the mass-weighted velocity of the gas contained in 
this sphere. 



3.1.2 Pop II star formation 

W e treat Pop I I star formation with the same prescription 
as IWise fc Cenl (|2009h . which is similar to the Pop III pre- 
scription but without the minimum H2 fraction requirement. 
This is removed because the metal-enriched gas can effi- 
ciently co ol even in the presence of a strong UV radiatio n 
field (e.g. ISafranek-Shrader. Bromm fc Milosavlievic1 f2010). 
To ensure the volume is cooling, we restrict star formation 
to gas with temperatures T < 1000 K. Unlike Pop III star 
particles that represent individual stars, Pop II star particles 
represent a star cluster of some total mass and an assumed 
Salpeter IMF. 

Once a cell meets these criteria, the prescription 
searches outward with increasing radius for the boundary of 
the molecular cloud, centred on the most dense cell. Here, 



we define a molecular cloud as a sphere with a dynami- 
cal time £dyn = 3 Myr (corresponding to an average den- 
sity pmc — 500/i cm -3 ) and a radius Rmc, where /1 is 
the mean molecular weight. This sphere typically has a ra- 
dius of 6 pc (increasing as M^[q) for the smallest molecular 
clouds and 2,500 cells. Once the sphere radius is found, a 
fraction c* = 0.07/ co id of the cold gas (T < 10 3 K) in- 
side the sphere is converted into a star particle with mass 
m* = c*(47t/3)pmc#mc- This treatment of cold gas ac- 
cretion is similar to a star formation model with a multi- 
phase subgrid model (|Springel fc Hernquistl I2003T I but is 
employed in a simulation that can resolve the multi-phase 
interstellar stellar medium (ISM). After the star particle 
is created, we replace the sphere with a uniform density 
Pmc = (1 — c*)/(Gtd yn ) and temperature T = 10 4 K, which 
approximates the initial stages of an H 11 region. 



3.2 Stellar feedback 

The radiation field is evolved with adapt ive ray tracing 
i|Abel fc Wanderd |2002| ; IWise fc Abel 1201 j) that is based 
on the HEALPix framework ( Gorski et al.l I2005T ) and is 
coupled self-consistently to the hydrodynamics. Each star 
particle is a point source of hydrogen ionising radiation 
with the ionising luminosity equally split between 48 initial 
rays (HEALPix level 2). We use a mono-chromatic spec- 
trum for the radiation with the energy -B p h, the spectral 
shape weighted photon energy. For a cosmological simula- 
tion that focuses on galaxies, this simplification does not 
signifi cantly affect the o verall galactic dynamics (see Sec. 

6.3 in lWise fc AbellbOllT ). We do not consider He 1 and He 
11 ionising radiation. As the rays propagate from the source 
or into a high resolution AMR grid, they are adaptively split 
into 4 child rays, increasing the angular resolution of the so- 



lution, when the solid angle of the ray 9 — 4n/ '(12 x 4 ) is 
larger than 1/3 of the cell area, where L is the HEALPix 
level. We limit the ray splitting to a maximum HEALPix 
level of 13, or equivalently 8.05 x 10 s rays per source. We 
model the H2 dissociating radiation with an optically-thin, 
inverse square profile, centred on all Pop II and III star par- 
ticles. 



3.2.1 Population III stellar feedback 

We use mass-dependent hydrogen ionising and LW pho- 
ton luminositie s and lifetimes of the Pop III stars from 
ISchaererl (2002). We use a mass- independent photon energy 
E p h = 29.6 eV, appropriate for the near-constant 10 5 K 
surface temperatures of Pop III s tars. They die as Type II 
SNe if 11 ^ M^/M p, < 40 M i|Wooslev fc Weaver! Il995l ; 
iNomoto et al.ll2006l) and as PISNe if they are in the mass 
range 140-260 M (jHeger et ai]|2003h . where M* is the stel- 
lar mass. For normal Type II SNe between 11-20 Mq, we 
use an explosion energy of 10 5 1 erg and a linear fit to the 
metal ejecta mass calculated in lNomoto et al.l (|2006), 

Mz/Mq = 0.1077 + 0.3383 x (M*/M - 11) (10) 

or equivalently, the metal yield fraction 

y = 0.3383 - 3.614AC 1 . (11) 

We model the SNe of stars with 20 < M*/Mq < 40 as hy- 
pernova with the energies and ejecta masses (y ~ 0.15 — 0.2) 
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also taken from iNomoto et all linearly interpolating their 
results to M+. For PISNe, we use the explosion energy from 
iHeeer fc Wooslevl (|2002h . where we fit the following function 
to their models, 



Er- 



= 10 51 x [5.0 + 1.3O4(M Hc /M - 64)] erg, 



(12) 



where M Hc = (13/24) x (M* -20)M© is the helium core (and 
equivalently the metal ejecta) mass and M* is the stellar 
mass. If the stellar mass is outside of these ranges, then an 
inert, collisionless black hole (BH) particle is created. 

The blastwave is modelled by injecting the explosion 
energy and ejecta mass into a sphere of 10 pc, smoothed 



at its surface to improve numerical stability (jWise fc Abell 
2008). Because we resolve the blastwave relatively well with 
several cells across at its initialization, the thermal energy 
is converted into kin etic energy and a grees with the Sedov- 
Taylor solution (e.g. iGreif et al.ll2007t ). 

3.2.2 Population II stellar feedback 

Pop II star particles emit 6000 hydrogen ionising pho- 
tons per stellar baryon averaged over their lifetime and 
-Eph = 21.6 eV, ap propriate for a [Z/H] = —1.3 popula- 
tion (jSchaererl 12003). The star particles live for 20 Myr, the 
maximum lifetime of an OB star. These stars generate the 
majority of the ionising radiation and SNe feedback in stellar 
clusters, thus we ignore any feedback from lower mass stars. 
By considering a constant luminosity, we may be underesti- 
mating the impact of radiative feedback because, given an 
IMF, the total luminosity of a cluster is maximal at early 
times when thermal and radiation pressure forces are the 
greatest, i.e. the edge of H II region is still near the massive 
stars. 

For SN feedback, these stars generate 6.8 x 10 48 erg s _1 
Mq -1 after living for 4 Myr, which is injected into spheres 
with a radius of 10 pc. However, if the resolution of the grid 
is less than 10/3 pc, we distribute the energy into a 3 3 cubic 
region surrounding the star particle. Here, star particles also 
return ejected material with mass 

0.25 At Mi, /JW 
Am ej = - _ - M ^ , (4 Myr <t-t hilth <t ) (13) 

back to the grid at every timestep on the finest AMR level, 
to = 20 Myr denotes the lifetime of the star particle. This 
ejected gas has solar metallicity Z — 0.02, resulting in a 
total metal yield y = 0.005. However, we note that cal- 
culations which are calibrated to Milky Way observations 
jMadau et al.lll996l ) suggest y ~ 0.02, and thus we may be 
underestimating the Pop II chemical feedback. 

3.3 Simulation variations 

All three simulations start from identical 1 Mpc 256 3 initial 
conditions, have the same refinement criteria, and use the 
same star formation and feedback models. Starting with our 
reference model in Paper I, each following simulation builds 
upon the physics of the previously described one. 




15 20 
Redshift 

Figure 2. Time evolution of the Lyman- Werner background in 
uni ts of J21 = 10~ 21 e rg s~ 1 cm - 2 Hz -1 sr -1 , calculated with 
the I Wise fc Abell j2005fl model. See eq. ^ for its functional fit. 



in enzo ()Abel et al.l 1 19971: lAnninos et all 19971 ) and the 
H2 cooling rates from iGlover fc Abell (2008). We spatially 
distinguish metal enrichment from Pop II and Pop III stars, 
however they do not contribute to the radiative cooling rates 
and act as passive tracer fields for the SNe ejecta. We will 
refer to this simulation as the "Base" simulation. It required 
150,000 CPU hours on 512 compute cores to run to z = 8.1. 



3.3.2 Metal cooling 

In this simulation, we only add the effects of radiative cool- 
ing from fine-structure transitions in metals. We will re- 
fer to this simulation as the "MC" simulati on. We use the 
method of ISmith. Sigurdsson fc Abel (|2008l ) that is incor- 
porated into enzo v2.0. We calculate the cooling rates with 
Cloudy 07.02.01, which has a sophisticated chemical net- 
work that includes all atomic species and many molecular 
specie^], on a logarithmically spaced 4-dimensional grid of 
(i) density - 10" 6 to 10 6 cm" 3 , A = 0.25 dex, (ii) temper- 
ature - 10 to 10 s K, A = 0.1 dex, (iii) electron fraction - 
10" 6 to 1, A = 0.25 dex, and (iv) metallicity - 10" 6 to 1 
Zq , A = 1 dex, where A is the spacing between grid points. 
We use the CMB radiation spectrum as the input spectrum. 
Dust cooling becomes important at > 10 9 cm~ 3 , which our 
simulations are not designed to probe, and thus we do not 
consider cooling from dust grains. It required 200,000 CPU 
hours on 512 compute cores to complete. Compared to the 
Base simulation, the additional compute time comes from 
ray tracing through a larger ionised volume, i.e. the rays are 
not terminated when fully absorbed. 



3.3.3 Radiation pressure and soft UV background 

This simulation adds momentum transfer from ionising 
radiation to gas and an H2-dissociating radiation back- 



3.3.1 Reference model 

We use the nine-species (H I, H n, He I, He II, He 
III, e~, H2, H2 + , H~) non-equilibrium chemistry model 



2 Because the cooling rates from rotational and vibrational 
H2 transitions are already included in ENZO, thus we do not in- 
clude H2 cooling in the Cloudy rate table. 
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ground. We will refer to this simulation as the "RP" sim- 
ulation, and the runtime is similar to the Base simulation 
at 150,000 CPU hours. As the radiation field is calculated 
with enzo+moray, the momentum of absorbed ionising ra- 
diation 



dPr P = r 



(14) 



is transferred to the absorbing medium, where dP is the 
number of photons absorbed, E p ^ is the photon energy, and 
f is the normal direction of the radiation. This additional 
momentum further accelerates the gas by 

at p Vcell 

per unit mass, where p is the gas density, dt is the radiative 
transfer timestep, and V ce ii is the cell volume. Currently we 
do not consider momentum transfer by radiation trapping 
when UV light is absorbed by dust grains and re-emitted in 
the IR and absorbed many times. In principle, this would 
increase the da ap by a factor of /trap, which we discuss in 
Section 15.21 The acceleration da ap is then added to the ac- 
celeration field in the calculation, for example, from gravity, 
in an operator split fashion. 

A time-dependent Lyman- Werner (LW) optically thin 
radiation background is utilised in this simulation on which 
the LW radiation from point sources are added. LW flux 
from point sources dominates over the background within 
3 kpc, dissociating H2 and thus delaying Pop III star for- 
mation in nearby h a loes llMachacek. Bryan fc A bel 2001; 
IWise fc Abel l2007bl : lO'Shea fc Normanl |2008h . To calcu- 
late its intensity as a function of time , we use the semi- 
analytical model of I Wise fc Abell |2005), updated with the 
7-year WMAP parameters and optical depth to Thomson 
scattering. The LW radiation intensity is plotted in Fig. 
In this semi-analytical model, we use a Pop III stellar mass 
of 100 Mq, st ar formation efficie ncy of 0.005, and escape 
fraction of 0.2 (Wise & Ccn 2009). The intensity decreases 
after z ~ 14 because Pop II star formation becomes dom- 
inant, which produce less LW specific luminosity than Pop 
III stars that have surface temperatures of T = 10 5 K. For 
computational convenience, we fit the background evolution 
with the function 



logio J2i(z) 



A + Bz + Cz 2 + Dz 3 + Ez 4 , 



(16) 



where (A, B, C, D, E) = (-2.356,0.4562,-0.02680,5.882 x 
10~ 4 , —5.056 x 10 -6 ) and J21 is the specific intensity in units 
of 10 -21 erg s _1 cm -2 Hz -1 sr -1 . This fits the model data 
within 1 per cent in 6 ^ z ^ 30 and is consistent with J21 
values in lTrenti fc Stiavellil (|2009h . 



4 RESULTS 

We concentrate on the most massive halo in the simulation 
to demonstrate the effects of momentum transfer from ionis- 
ing radiation. For reasons discussed later, the halo in the MC 
simulation overcools and forms stars at high specific SFR, 
which ionises the entire simulation volume at z = 8.1. Solv- 
ing the radiative transfer equation with enzo+moray in 
this simulation becomes computationally unfeasible in the 
optically thin limit, so we stop the MC simulation here and 
compare the halo in the three simulations that this point. 



At redshift 8.1, the halo has a DM mass of 2.0 x 10 s M©, 
corresponding to a virial temperature 1.7 x 10 4 K and cir- 
cular velocity of 21 km s _1 . This halo experienced a ma- 
jor merger at z ~ 10. It is undergoing a major merger at 
2 = 8.1 that increases the virial temperature to > 10 4 K at 
which point it can efficiently form stars. The gas fraction of 
the halo f g = Mg as /M to t has recovered to 0.13, previously 
being depleted to under 5 per cent in the Pop III hosting 
progenitors that had its gas expelled by stellar and super- 
nova feedback. Here M to t is the total mass (including DM) 
of the halo. At this time, the Base, MC, and RP simulations 
have a stellar mass of 2.5x 10 5 , 6.7 x 10 5 , and 4.5 x 10 5 M , 
respectively. By redshift 7, its total mass is 1.0 x 10 9 M Q , 
containing 2.1 x 1O 6 M0 of stars in the Base simulation. We 
refer the reader to Paper I for more details at z = 7. 



4.1 Comparison of star formation histories 

The top panel in Figure [3] shows the evolution of the ionised 
volume fraction of all three simulations. The second, third, 
and fourth panels show the star formation rate (SFR). The 
specific SFR (sSFR), the gas accretion rate M gas , and the 
ratio of SFR to M gas are shown in the bottom three pan- 
els. At redshift 12, both MC and RP simulations produce 
metal-enriched stars after the halo r e aches the filtering m ass 
( Gnedin fc Huilll998l ; ICnedinl [20051 ; IWise fc Abeil 120081 ) of 
~ 1O 7 M0. The burst in the Base simulation is slightly de- 
layed and weaker by an order of magnitude because of the 
lack of metal cooling. In all three simulations, about 15 per 
cent of the accreted gas is converted into stars. Radiative 
and supernova feedback in the stronger bursts of the MC 
and RP runs drive a net outflow from the halo, seen in neg- 
ative values of Af gas . 

The halo continues to accrete from numerous minor 
mergers of ~ 10 6 Mq haloes until z — 8.1. One notable 
difference in the RP simulation is the enhanced gas ac- 
cretion between redshift 9 and 10. Here, the UV back- 
ground suppresses Pop III star formation in minihaloes 
by increasing the minimum mass to cool and condense 
by H 2 alone dMachacek. Bryan ^fc Abell l200ll ; IWise fc Abel 
l2007bl ; lO'Shea fc Normanl I2008T ). This results in more gas- 
rich mergers, thus increasing the gas accretion rate. The halo 
continues this constant but slow mass buildup until a major 
merger at z = 8.1. Before the major merger, the SFRs are 
similar between 10 -2 and 10 -3 M yr -1 in the simulations 
with the RP simulation sustaining a more consistent SFR 
mainly because of the additional gas accretion whereas the 
other simulations have more a bursty behavior. The galaxies 
convert between 0.1 and 1 per cent of the accreted gas into 
stars during this quiescent period before it reaches a virial 
temperature of 10 4 K. 

The most significant difference between the simulations 
is the star formation when the halo reaches a virial temper- 
ature of 10 4 K through a major merger. The Base simula- 
tion has a SFR = 3.7 x 10~ 3 Mq yr -1 and a sSFR = 15 
Gyr -1 . In the MC simulation, the additional radiative cool- 
ing provided by the metal lines prompts a SFR = 9.9 x 10 -2 
Mq yr -1 and a sSFR = 150 Gyr -1 , an order of magnitude 
higher than the Base simulation. The radiative and super- 
novae feedback from these stars drives a net gas outflow of 
-0.2 Mq yr -1 from the halo, as seen in the fourth panel of 
Figure[31 but cannot abate the star formation alone. Momen- 
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Figure 4. A visual comparison of the most massive (2.3x 10 8 Mq) halo at z = 8.1 in simulations with the base physics model (left column), 
including radiative cooling from metal species (middle column), and including radiation pressure and an H2-dissociating background 
(right column). The different rows contain density- weighted projections of gas density (top), temperature (middle), and total metallicity 
(bottom). The field of view is two virial radii (3.4 kpc). Comparing the Base and Metal Cooling simulations, the metal-rich ejecta cannot 
expand before all of its thermal energy is radiated away, primarily through the additional metal cooling. The additional impetus from 
radiation pressure alleviates this artifact and can regulate star formation by driving gas outflows and further entrainment with the 
surrounding medium. 



turn input from massive stars in the RP simulation reduces 
the SFR by a factor of five to 1.8 x 10" 2 M yr -1 and the 
sSFR to 37 Gyr -1 . On average, approximately 3 per cent 
of the accreted gas forms stars in the latest star formation 
period, and the variations are caused by stellar feedback. At 
2 = 8.1, the enhanced star formation in the MC simulation 
ionises 96 percent of the simulation volume, which is in stark 
comparison with the Base and RP simulations with ionised 
fractions of 42 and 36 percent, respectively. The RP ionised 
fraction is slightly lower because of the suppression of Pop 
III SF from the LW background. 

Observations of star-forming galaxies show a tight re- 
lation between sSFR and stellar mass, and we can compare 
our results to determine our most physical model. Massive 
(M* ^ 1O 1O M ) galaxies at z ^ 2 have sSFR a M~ 0A , 
flattening to an upper thres hold of 2.5 Gyr - in low-mass 
galaxi es (|Karim et al.l 120111 ). At z = 4 - 6, IStark et ail 
(2009) found a strong sSFR-M* correlation with a simi- 



larly large scatter (see Figure 1 in iKhochfar fc Silk! 1201 ll ). 
In that work, the average SFR is ~ 3 Mq yr -1 in sys- 
tems between M* = 10 8 and 10 11 Mq, or equivalently, an 
average of 30 Gyr -1 in stellar systems of 1O 8 M . Galax- 
ies at z > 7 also show a large scatter in sSFR, ranging 
from 1 to 100 Gyr -1 for galaxies with stellar masses be- 
tween 10 7 and 1 9 Mq, trending toward higher sSFRs in 
smaller systems dSchaerer fc de Barrosl |2010| ; iMcLure et al.l 
l201ll IKhochfar fc SiUdl201ll ). Do high-redshift galaxies have 
more gas available for star formation to sustain sSFR > 2.5 
Gyr -1 , as suggested by the z > 7 data? Does this relation 
continue to the very smallest galaxies at high-redshift? If 
we assume so, the MC simulation produces the best match 
with the other models falling below the relation. This goes 
against intuition in that a more realistic model, in our case, 
adding radiation pressure and a LW background, produces 
more realistic galaxies. To further understand this behav- 
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Figure 5. Mass- weighted average metallicities as a function of density and temperature in the most massive halo at z = 8.1 for the 
Base (left), metal cooling (middle), and radiation pressure (right) simulations. When metal cooling is considered, metal-rich gas above 
3 X 10 -23 g cm -3 can cool to the CMB temperature. The majority of gas with > 10 -2 Zq is restricted to SN remnants (hot ISM phase 
up to T = 10 7 K) and cold and dense gas. The metals are not mixed into the lower density ISM. Radiation pressure (right panel) 
provides an additional impetus to mix the metal-rich ejecta with the entire ISM, resulting in nearly all of phase-space being enriched to 
an average of ~ 10 — 2 Zq. 



ior, we must investigate the gaseous properties and stellar 
population of the galaxy. 

4.2 Overcooling and artificially enhanced star 
formation 

Here we concentrate on the causes behind the enhanced star 
formation in the metal cooling simulation when compared to 
the Base simulation. Figure [4] shows density- weighted pro- 
jections of gas density, temperature, and total metallicity 
for the three simulations. The gas distributions are irregular 
and have no organised disc rotation in any of the runs. There 
are three DM filaments t hat accrete onto the halo, which is 
typic al of rare peaks fe.g-. lDekel et al' I l2009l ;l Danovich et al.l 
l201ll ). Only one filament has retained its gas, which is lo- 
cated to the left of the halo in Figure [3] The other filaments 
have been photo-evaporated and Jeans-smoothed by stellar 
radiation. As previously mentioned, the MC simulation is 
nearly fully ionised with an average intergalactic medium 
(IGM) temperature of ~ 3 x 10 4 K, clearly seen in the tem- 
perature projections. The additional pressure provided by 
the heated IGM compresses the filament in comparison with 
the Base and RP simulations. 

The spatial metallicity distributions are the most rele- 
vant in understanding the SFR deviations between the simu- 
lations. In the Base simulation, the metal-rich outflows with 
[Z/H] between -1 and -2 can be seen expanding out to a 
radius of ~ 1 kpc. The MC simulation produces a factor 
of three more stars, but the metal-rich ejecta are restricted 
in a small volume with a radius of ~ 0.5 kpc with a mean 
metallicity of [Z/H] = -0.5. These metal-rich regions are 
surrounded by a metal-poor interstellar medium (ISM) with 
-4 < [Z/H] < -2. 

Figure [S] shows the metallicity as a function of density 



and temperature, and it further illustrates the differences 
that metal cooling induces. Without metal cooling in the 
Base simulation, the gas can only cool to ~ 1000 K through 
H2 cooling in the gas phase. The metals are mixed to an 
average of 3 x 10~ 2 Z in gas with p > 3 x 10~ 26 cm~ 3 . The 
nearly solar metallicity isobars reaching 10 7 K are Type II 
SNe remnants. When metal cooling is considered, the rel- 
atively metal-rich (10 _1 Zq) material is mainly confined to 
the cold phase of the ISM above 1 cm -3 and below 1000 K 
down to the CMB temperature. This cold, metal-rich gas 
mixes little with the surrounding diffuse ISM that has a 
mean metallicity of 10 -3 Zq. 

Because stars form in cold and dense gas, the stellar 
metallicities reflect the gas metallicities in the cold phase. 
They can be a plausibility check as they can be indirectly 
compared to observed metallicities in local dwarf galaxies, 
which have a clear luminosity-metallicity relationship that 
increases from [Fe /H] = -2.5 to -1.5 over a luminosity range 
log(L/L Q ) = 3-7 (|Kirbv et alj|201ll ). 

Figure|5]compares the metallicity distribution functions 
(MDFs) in the most massive galaxy at z = 8.1 in each 
simulation. The Base simulation has a bimodal distribution 
within the range —2.5 < [Z/H] < —1.0. The lower metal- 
licity stars were formed before the halo reaches T v i r = 10 4 
K. Afterward, the stars quickly enrich the halo to [Z/H] ~ 
-1.5, which forms stars at similar metallicities, creating a 
bimodal MDF. This process is described in more detail in 
Paper I. A similar enhancement in star formation occurs in 
the MC simulation at T v i r = 10 4 K but to a greater ex- 
tent because of higher cooling rates and thus star formation 
rates. Stars primarily form in the inner 500 pc, where the 
gas is enriched to [Z/H] ~ -0.5 because the metal-rich ejecta 
is confined within these inner regions. 

These metallicities are a factor of 30 above the 
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Figure 3. The top panel shows the ionised volume fraction, while 
the next three panels show the star formation history for the 
reference, metal cooling, and radiation pressure simulations. The 
fifth panel displays the specific star formation rate for all three 
simulations. The sixth panel shows the gas accretion rate in the 
halo, and the bottom panel shows the ratio of star formation and 
gas accretion, which is effectively the star formation efficiency, 
fncluding metal cooling results in overcooling in the halo interior, 
efficiently producing stars at a rate of ~ 0.1 yr — 1 . This single 
burst ionises the entire simulation box. Radiation pressure from 
ionising radiation regulates star formation to ~ 0.01 Mq yr — 1 . 
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Figure 6. Normalised stellar metallicity distribution functions 
of the most massive halo for the Base (solid lines), metal cool- 
ing (dashed), and radiation pressure (dash-dotted) simulations at 
z = 8.1. Overcooling in the metal cooling simulation creates a 
stellar population that is nearly solar metallicity, whereas the ad- 
ditional mixing with radiation pressure shifts the average back to 
~ 10~ 2 Zq, agreeing with the mass-metallicity relation of galax- 
ies. 

luminosity-metallicity relationship and suggest that our MC 
model produces galaxies that are too luminous and metal- 
rich. The extra radiative cooling from metal species allows 
for the super-solar SNe ejecta to lose most of its thermal en- 
ergy before it launches a high-velocity outflow, and most of 
the ejecta are trapped within the galaxy. Then this metal- 
rich material cools again and condenses into star forming 
clouds. This process repeats several times, creating a run- 
away effect in overproducing stars, known as the well-known 
and studied numerical "overcooling problem" that produces 
gal axies that are too centrally concen trated (first reported 
bv lKatz. Weinberg fc Hernquist|[l99rj ). This was solved by 
introducing SNe feedback and conducting the simulations 
at higher resolutions; however, our calculations show the 
same problem at parsec-scale spatial resolution and despite 
resolved SNe feedback. 



4.3 Effects of radiation pressure 

We find that momentum transfer from ionising radiation to 
the absorbing gas alleviates the overcooling problem seen in 
the MC simulation. Recall that at z — 8.1, the SFR is re- 
duced by a factor of five compared to the MC simulation. To 
understand this difference, we focus on its effect on internal 
gas dynamics and the chemo-thermal state of the gas, the in- 
creased ejecta dispersion from their origin, and the ensuing 
star formation. 



4-3. 1 Gas dynamics 

Unlike the thermal energy produced by SN feedback, mo- 
mentum must be conserved and is not lost to radiative pro- 
cesses. The momentum transferred to the absorbing gas gen- 
erates a medium in the central 100 pc that is both tur- 
bulent and expanding radially, where the young stars and 
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their H II regions are clustered. As the D-type shocks ex- 
pand with the ionisation fronts, momenta from the absorbed 
photons are transferred to the optically thick shock. Their 
cumulative action produces an average outflow in the inner 
100 pc, which is visible in the radial profile of the mass- 
averaged radial velocity, shown in the top panel of Figure [7| 
The radial velocity is calculated after subtracting the mass- 
averaged bulk velocity of the halo. In the central 15 pc, the 
mass-averaged radial velocity is 20 km s _1 , and it gradu- 
ally decreases to 2 km s _1 at 40 pc. This weak expansion of 
the inner dense gas extends out to 100 pc. In contrast, the 
Base and MC simulations both exhibit weak inflows up to 
5 km s~\ where SN feedback is not enough to support the 
dense central gas. 

Within the inner 100 pc, there are 49 young (< 20 Myr) 
star particles with a total mass of 4.6 x 1O 4 M . The interac- 
tion between the expanding shocks associated with their H 
n regions and SN blastwaves creates and sustains this tur- 
bulent medium. To measure the turbulent motions, we cal- 
culate the three-dimensional rnis velocity u rm8 as a function 
of radius. It is computed with respect to the mass-averaged 
velocity of each thick shell. On average, it is ~ 20 km s _1 
and varies little with radius, indicating that turbulence is 
widespread throughout the halo. The turbulence is mildly 
supersonic at r > 100 pc with turbulent Mach numbers 
Aft ur b = (v T ms) /(c s ) decreasing with radius from two at the 
virial radius, where c s is the sound speed and the angled 
brackets denote a mass-weighted average in the shell. This 
is nearly the circular velocity of the halo, V c = 21 km s" 1 . In 
the inner 50 pc, v rmB ~ 20 km s _1 , and it slowly decreases to 
15 km s _1 at 130 pc and increases to 20 km s~ at the virial 
radius. Outside 300 pc, radiation pressure has little effect 
on turbul ent motions, where i t is primarily dri ven by viri- 
alization l|Wise fc Abelll2007al ; iGreif et al.l [20081. The Base 
simulation drops to 10 km s _1 within r = 100 pc as the 
gas cools and condenses. The loss of gaseous kinetic energy 
is even more apparent with metal cooling, lowering u rms to 
5 km s _1 in the inner 50 pc. This is consistent with an over- 
cooling core that is not pressure supported by either thermal 
or turbulent energy. 

Interestingly, the turbulent Mach numbers, shown in the 
bottom panel of Figure[7] are almost identical in the MC and 
RP models in the inner 100 pc. This occurs even though the 
addition of momentum feedback creates a « rms four times 
larger than the MC simulation. This behavior is caused two 
complimentary processes. First, radiation pressure creates 
a more diffuse medium in which the SN explosions occur, 
which results in a lower cooling rate and keeps the tem- 
peratures and, equivalently, sound speeds high. Second, the 
additional driving from radiation pressure and blastwaves 
that have not stalled in the ISM keep the turbulent veloci- 
ties nearly equal to the halo circular velocity. These physical 
processes lead to an unchanged turbulent Mach in both star 
forming environments. 

The velocity dispersion is a good measure of the gen- 
eral strength of turbulent motions, but the vorticity u> better 
characterises the local strength of turbulence. In Figure [8j 
we plot two-dimensional histograms of vorticity (top panels) 
as a function of radius and its growth (bottom panels), which 
is the right hand side of Equation For the RP simula- 
tion, we include the baroclinic term from radiation pressure 
in our analysis, (a rp x Vp) /p. To visualise the fluid flow, it is 




Radius [pc] 



Figure 7. Top panel: Mass-weighted radial profiles of gas radial 
velocity of the most massive halo, centred on the halo's centre 
of mass, from the reference (solid), metal cooling (dashed), and 
radiation pressure (dash-dotted) simulations. The Base and MC 
runs exhibit a zero to small infall velocity within 300 pc, whereas 
adding radiation pressure results in an average gas outflow in the 
inner 100 pc. Middle panel: Three-dimensional rms velocity as a 
function of radius, illustrating the turbulent medium created by 
momentum transfer from ionising radiation. Turbulence motions 
support the object from further collapsing and assist to mix the 
metal-rich ejecta with the metal-poor ISM. Bottom panel: Radial 
profile of the turbulent Mach number Afturb' The MC and RP 
simulations exhibit nearly identical Mt uf b values in the inner 100 
pc, showing that the turbulence decays in the MC run as the gas 
rapidly cools. 



helpful to consider the local rotational period is 47r/|u>|. Gas 
with the highest vorticities exist at 100 < r/pc < 500. This 
is primarily driven by virial shocks as matter is accreted 
from the IGM and filaments. This region is approximately 
where the density increases by an order of magnitude in the 
density projections of Figure|4l The MC and RP simulations 
have similar average vorticities until r < 50 pc, where the 
RP simulation exhibits an increase of 10, whereas the MC 
simulation stabilises at 2 x 10~ 14 s _1 . Vorticity generation is 
significantly higher by a factor of 100 within this inner region 
with radiation pressure. We find that non-barotropic flows 
from thermal pressure are primarily responsible for sustain- 
ing these turbulent motions. The contribution from vortic- 
ity stretching and radiation pressure provides ~10 _1 and 
~10~ 5 of Dtj/Dt, respectively. The very weak driving from 
radiation pressure does not necessarily indicate that it can 
be neglected. Because it adds momentum to the shells, the 
shocks during shell collisions are stronger and produce more 
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Figure 8. Top row: Mass-weighted histograms of the magnitude 
of vorticity for the Base, MC, and RP runs from left to right. The 
mass- weighted radial profile of vorticity is shown by the solid line. 
Bottom row: Total time derivative of vorticity, i.e. right hand side 
of Equation JgJ. Within the inner ~70 pc, vorticity is generated 
at a higher rate in the RP simulations than the Base and MC 
cases, leading to a more turbulent medium that is also depicted 
in the rms velocity (Figure [7J. 



vorticity through the thermal baroclinic term. The weak 
contribution from radiation pressure is straightforward to 
understand. It is strongest when the ionisation front radius 
is small, and the H n is close to spherical. In the spherically 
symmetric case, radiation will always be absorbed parallel 
to the density gradient as the optical depth increases, and 
thus a rp x Vp is small. In a clumpy medium, this term will 
be larger when radiation grazes the surface of dense clumps. 
It should be noted that by neglecting scattering we could 
be underestimating its magnitude, and our calculations are 
conservative in quantifying the impact of radiation pressure 
in the ISM of dwarf galaxies. 



4-3.2 Metal mixing 

The gas metallicity shown in Figure [3] show that the SN 
ejecta has expanded beyond the virial radius with outflows 
with velocities up to 250 km s" 1 , enriching most of the gas 
to [Z/H] ~ -2. This behaviour is further illustrated in Figure 
where the SN ejecta thoroughly mixes with the ISM at all 
density-temperature pairs in the cold and warm phases. This 
greatly differs from the lack of mixing in the MC simulation 
that creates a distinct separation between dense metal-rich 
gas and diffuse metal-poor gas. 

Metal mixing is enhanced by a sustained turbulent 
medium, resulting in a MDF (Fig. [5| with a mean metal- 
licity of [Z/H] = —2.1 and a standard deviation of 0.2 
dex. This is in excellent agreement with the luminosity- 
metallicity relationship f or dwarf galaxies fo r a stellar system 
of A'L = 4.5 x 10 5 M© (|Kirbv et alj|201ll ). This agreement 
suggests that radiation pressure plays an important role in 
regulating star formation in high-redshift dwarf galaxies. 



4-3.3 Star formation regulation 

We have demonstrated that momentum input from massive 
stars prevents the overcooling and over-enrichment of the 
gas in galaxy formation simulations. The additional momen- 
tum from radiation both disperses the metal-rich SN ejecta 
to large radii and provides turbulent support to the system. 
The first effect reduces the metallicity in the dense gas. This 
reduces the radiative cooling from metal species and the 
amount of cold, dense gas available for star formation. The 
latter effect prevents the gaseous core from catastrophically 
collapsing and forming stars at a high efficiency. 

Our comparison between the MC and RP simulations 
shows that star formation is self-regulated with radiation 
pressure playing a vital role in early galaxy formation. To 
reiterate, the RP simulation creates a galaxy with a sSFR = 
37 Gyr -1 and a stellar population whose MDF ({Z/Zq) = 
—2.1) is in excellent agreement with observations of local 
dwarf galaxies. 



4.4 Radiation pressure characteristics 

Here we describe some of the properties of radiation pressure 
within our largest galaxy. Figure [5] shows two-dimensional 
slices of acceleration from momentum transfer a lp (Equation 
IT5]) . gas density, temperature, and total metallicity through 
the centre of mass of the galaxy. The acceleration a rp is 
scaled with respect to the gravitational acceleration at the 
virial radius, a gr av = GM viv /rl ir = 7.0 x 10 -10 cms -2 . The 
time to accelerate the gas above the escape velocity v csc — 
(2GM vil /r vil ) 1/2 is 



tcsc = fosc/a g rav = [50fi m (z) R 1 (z)] 1/2 



(17) 



which is independent of the halo. Here Q m (z) and H(z) 
are the matter density and Hubble constant at redshift z. 
In a matter dominated universe (z 3> 1) for the adopted 
ACDM cosmology, this timescale is 120[(l + z)/10]~ 3/2 Myr, 
which is possible with sustained star formation in a galaxy. 
However, this should be only be used to gauge the instan- 
taneous impact of the momentum transfer because the ac- 
celeration will decrease as the shell travels farther from the 
source. 

The density plot shows that the ISM is fractal and 
porous in nature with the H n regions, creating dense shells 
on their surfaces. The acceleration a lp exceeds 10 per cent 
of flgrav in most of the H n region except for the interiors 
of the SN remnants, where the neutral hydrogen fraction 
is nearly zero. It exceeds a gra v in the dense, neutral shells 
surrounding the H II regions. Figure [TOl shows the probabil- 
ity and cumulative distribution function of a rp /a gI . a v inside 
of this halo. There is a plateau in the distribution between 
2 x 10~ 3 and 0.3 at 2 x 10 6 Mq per dex. Radiation pressure 
is accelerating about 4 per cent (~ 3 x 10 s Mq) of the halo 
gas greater than a gra v 

The slices of temperature and total metallicity in Figure 
[9]show the hot phase of the ISM, which have T > 10 6 K and 
approximately solar metallicity. Our stellar feedback pre- 
scription injects SN thermal energy and metals after 4 Myr 
of main sequence radiative feedback. Consequently, the SN 
remnants are initially contained within the rarefied medium 
of the interior of their associated H II regions and then ex- 
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Figure 9. Top panel: Slices of instantaneous acceleration from 
momentum transfer from absorbed ionising radiation, which has 
been scaled to the total gravitational acceleration of the halo, 
Hgrav = GMvix/r^- . Larger values indicate a greater outward 
acceleration from the responsible star particle. The arrows denote 
the direction and magnitude of the acceleration field. The field of 
view is 1 kpc and is cut through the halo's centre of mass. Bottom 
panels: The same field of view for gas density (left), temperature 
(middle), and metallicity (right). 



pand into the ISM, and in some cases, they drive metal-rich 
outflows into the IGM. 



5 DISCUSSION 

Stellar feedback comes in various forms. First, photo- 
ionisation and photo-heating from massive stars create H 
II regions that provide pressure support to the warm phase 
of the ISM. An extreme case of photo-heating feedback 
exists in haloes with M < 1O 7 M0 that host Pop III 
stars, where the ~30 km s _1 shock launched by the H 
ii reg ion can expel over 90 per cent of the gas from the 
halo dWhalen. Abel fc Normanll2004 iKitavama et al.ll2004 
I Abel, Wise fc Bryan! I2007T T Second, SNe feedback injects 
thermal energy and heavy elements into the ISM, creat- 
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Figure 10. Mass-weighted probability density distribution (solid 
line) of the instantaneous acceleration from momentum transfer 
from absorbed ionised radiation, which has been scaled to a grav = 
GM v - lr /r 2 . = 7.0 X 10~ 10 cm s -2 , inside the virial radius of the 
most massive galaxy. The inverse cumulative distribution is shown 
by the dashed line. Four per cent (~ 3 X 10 5 Mq) of the gas is 
experiencing an acceleration greater than the total gravitational 
force of the halo. 



ing the hot bubbles an d blastwaves. Th ey can self-enrich 
their environment (e .g. Wise ct al. 20 l jf), d i srupt its star- 
form ing cloud (e.g. iDraine fc Woods! Il99ll ; IWhalen et alj 
l2Q08f), create realistic bulgeless d isc galaxies, and pos- 
sibly launch larg e-scale outflows l)Governato et all l2010l : 
iBrook et al.ll2012h . Both processes together can play an im- 
portant role in regulating star formation and altering gas 
and stellar metallicities, galactic morphologies, gas fractions, 
and the structure of the ISM. 

Another feedback process is the momentum transfer 
from ionising radiation from massive stars, i.e. radiation 
pressure. It has been suggested that radiation pressure can 
regulate s tar format i on an d launch galactic outflows in large 
galax ies (lHaehneltl 19951; Murray. Quataert fc Thompson! 



2005: 



Murray. Menard fc Thompson! 



Hop kins. Quataert fc Murray! l201ll ). In this pap 



2011; 



have investigated the role of momentum deposition from 
ionising radiation in the formation of high-redshift dwarf 
galaxies. We have found that the addition of radiation pres- 
sure regulates star formation through additional turbulent 
support from shell collisions and allowing for SNe to drive 
large-scale outflows from the galaxy, ejecting metal-rich gas 
in the process. These effects alleviate the galaxy overcooling 
problem seen in various simulations of galaxy formation. 

By including radiation pressure, our simulation creates 
realistic high-redshift dwarf galaxies that matches the local 
luminosity-metallicity relation. Although this is not a direct 
comparison between high-redshift and local dwarf galaxies, 
the metallicities of the stellar population is a good gauge of 
a reasonable ISM chemo-thermal state. If we assume that a 
fract ion of high-redshift galax i es survive until the pr esent- 
day (jGnedin fc Kravtsovll20o5 iBovill fc Ricottill201ll ). then 
the bounds of MDFs should not exceed the ones observed 
in local dwarf galaxies, which typically have most of the ir 
stars within the range -3 < [Fe/H] < -1 dKirbv et alJl201.ll ). 
For instance in our MC simulation, the stars have a median 
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metallicity of ~0.3 Zq, which will still exist in its MDF if we 
were to passively evolve them to z = 0. On the other hand, 
this M = 2 x 10 s Mq galaxy is most likely to be incorporated 
into a galactic bulge, b ased on its halo mass accretion history 
l|Wechsler et al.lf2002f) . However, this result of overcooling is 
generic to galaxies reaching T = 10 4 K and should occur in 
dwarf galaxies that form later. Another possible way to avoid 
this restriction is that these stars form from a top-heavy 
IMF. If all of the stars form with M > 0.8M©, then they 
will not survive until the present-day, and the MDFs from 
high-redshift dwarfs cannot be constrained by the observed 
MDFs in local dwarf galaxies. 



5.1 Feedback prescriptions 

Implementations of stellar feedback vary widely between dif- 
ferent works. Here we discuss these variations, and their re- 
spective advantages and disadvantages. As galaxy formation 
simulations reach higher resolution, indi vidual star form- 
ing clo uds are being resolved. For example, iGovernato et alj 
(|2010T) found that a bulgeless disc galaxy forms instead of 
a bulge-dominated galaxy when the star formation density 
threshold was increased by a factor of 1000 to 100 cm -3 . 
This occurred because the simulation allowed the gas to 
cool and clump into smaller associations instead of a large 
gas reservoir at the galaxy centre. When simulations resolve 
these clouds, more small-scale physical processes need to be 
accurately modelled, and the less previous "subgrid" models 
are applicable because of their original assumptions of star 
formation on the galactic scale. 

As discussed before, there are four main components 
to stellar feedback: photo-ionisation, photo-heating, momen- 
tum transfer from radiation, and SN explosions, in addi- 
tion to stellar winds and planetary nebulae. The tempo- 
ral separation of the radiative and SNe feedback is crucial. 
In the case where both mechanisms have the same time- 
dependency, the injection of thermal energy and heavy ele- 
ments are initially dumped into a dense and cold medium, 
unaffected by stellar radiation. If the SNe feedback is de- 
layed for some amount of time, in our case, the minimum 
lifetime of an O-star of 4 Myr, then an H II region is created 
surrounding the star(s), which lowers the ambient density 
and increases the temperature to ~10 4 K. An optically-thick 
shell forms around the H II region, and radiation pressure 
accelerates the shell in a ddition to pressure forces from the 
H II region (|Spitzerll978l ). After a few Myr, the massive stars 
start to explode in this warm and diffuse medium, instead 
of a cold and dense medium. Because cooling rates are oc p 2 , 
the SN remnants are less likely to overcool in an implemen- 
tation that delays SNe feedback from radiative feedback. 

There have been a few types of prescriptions of SNe and 
momentum injection to create galactic winds and to regulate 
star formation. For the most part, they have an unphysical 
basis but replicate the observed properties of outflows and 
cannot study the origin of the winds. 

In galaxy simulations, it has been historically difficult 
to resolve the multi-phase ISM. In such a simulation, SNe 
thermal energy is injected into a large region, usually at its 
resolution limit. Per timestep, the total thermal energy is 
only sufficient to heat the gas to T ~ 10 5 — 10 6 K, where the 
cooling rates are maximal. If the blastwave were resolved 
to parsec-scales, then the temperature in the Sedov- Taylor 



solution should be T ~ 10 8 K, where the cooling rates are 
orders of magnitude lower. One approach is to neglect radia- 
tive cooling within the blastwave, allowing i t to adiabatically 
expand (e.g. Thacker fe Couchmanl |2000| ; IGovernato et al] 
l2007l;lGuedes et al.ll20lC Using this method results in real- 
istic late-type galaxy morphologies, but it can deviate from 
the observed luminosity functions because they do not sup- 
press star formation enough, particularly in low mass galax- 
ies, and do not drive win ds at sufficiently high velocities and 
mass-loading factors (e.g. lGuo et al.l201Ch . We note that an- 
other approach is to sporadically inject thermal energy only 
after enough SNe thermal feedback has accumulated to heat 
the gas to T > 10 7 K, where the cooling times are high. Here 
eject a does not overcool and follow s the Sedov- Taylor solu- 
tion (|Dalla Vecchia fc Schavell2012h. 

Several group s dSpringel fc Hernquistl 120031; 
Qppenheimer fc Pavel 120081 ; Dall a Vecchia fc Schavel 



20081 ; ISales et alj|2010h have modelled momentum transfer 
from radiation by adding momentum, so-called "kick- 
ing", to nearby particles to a velocity that is comparable 
to the circular velocity V c of the halo. Similar pre- 
scriptions exist for kinetic fee dback from SN explosions 
l|Dalla Vecchia fc Schavel l2008h . The dependence on the 
choice of this kick velocity, a constant value, a constant 
multiple of V c , or halo mass dependent velocity, has been 
shown to replicate the galaxy luminosity func tion and mass- 
metallicity relation (|Finlator fc Pavel |2008T ). The kicked 
particles are decoupled from the hydrodynamics until they 
escape the virial radius. The authors acknowledge and warn 
that this implementation is phenomenological, and it is suc- 
cessful in quantifying the im portance of momentum-driven 
winds in MW-like galaxies. iHopkins. Quataert fc Murray! 
(2011) used a more physical model for momentum loading 
that does not decouple the winds from the hydrodynamics 
and determines the momentum transfer directly from 
the stellar luminosity. The momentum transfer occurs 
within <10 pc, where they explore two schemes. Their 
stochastic model randomly kicks these nearby particles to 
the local escape velocity of the star-forming clump. Their 
continuous model adds the momentum Ap = (1 + r) p Td)L/c 
at each timestep, where rj p ~ 1 and can vary to account 
to uncertainties. They find that both models produce 
similar results. However, they do not solve the radiative 
transfer equation, so they must assume that the radiation 
is uniform within each H II region. Radiation pressure from 
a non-uniform distribution of massive stars may result in 
some cancellation of the pressure forces in the centre of 
the a ssociation, while still providing an ov erall outward 
force ([Socrates. Davis fc Ramirez- Rui j I 2008T ). In contrast, 
iKrumholz fc Thompson! (|2012T ) showed that the assumption 
of a uniform radiation field could lead to an overestimate 
the momentum transfer. As ionisation front instabilities 
grow, optically-thin gaps are created in the shock front, 
allowing for most of the flux to escape in these gaps and 
lowering the average momentum transfer per emitted 
photon. 

Ultimately the location where momentum transfer oc- 
curs depends where the radiation is absorbed. To model 
this requires solving the radiative t ransfer equation. None of 
these models with the exception of IKrumholz fc Thompson] 
utilises such a solver, but there are several simulation codes 
with coupled radiation transport. Momentum transfer from 
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Figure 11. Column density of metals with T < 1500 K, where 
dust would not be sublimated and can contribute to the opac- 
ity and enhance momentum transfer from stellar radiation. The 
figure is centred on the galaxy centre of mass. 



ionising radiation should be straightforward to implement 
in these codes because the absorption coefficient is already 
calculated and that absorbed photon momentum should 
be transferred to the absorbing medium, which is how 
enzo+moray simulates momentum loading. 



5.2 Radiation trapping from dust and Lyman a 

In this paper, we have neglected momentum transfer from 
dust absorption, which then re-radiates in the IR and ab- 
sorbed again many times. Thus our results represent a lower 
limit of star formation regulation from radiation pressure. 
The effective cross section of dust in the Milky Way is 



1(T 



cm per hydrogen atom in the near-IR (e.g. 
lDraindl2003h . or equivalently, an opacity k^mw — 5 cm 
g _1 in terms of total gas mass. It should be noted that 
Rosseland dust opacity is temperature dependent, scaling 
as T 2 below 150 K and can vary by factors of a few at 
150K < T < 1500K. Above 1500 K, even the most sta- 
ble dust species are sublimated and the Rosseland opacit y 
drops by a few orders of magnitude (|Semenov et al.| [2003). 
If we assume the dust extinction properties do not change 
with metallicity or redshift and scale the MW values with 
metallicity, the cross section per hydrogen atom scales with 
metallicity, giving «y = (Z/ZqJk^mw- The optical depth to 
dust absorption is 



Td = / p Kd dl — Kd,) 



(18) 



where Nz is the column density of metals, and /i and 
ran are the mean molecular wei ght and mass of a hydro - 
gen atom. We take Z = 0.023 (Grevesse fc Sauval|[l99Sh . 
Figure [11] shows the metal column density Nz of cold 
(T < 1500 K) gas in the inner 500 pc of the galaxy. Dust 
in gas T > 1500 K is destroyed and is not considered in 
the plot. Using Kd,MW 
range N z = (lO^.lO 1 



5 cm 2 g _1 , the column density 
2 corresponds to a optical 
). In this case, the dust ab- 



cm 

-2.7 



depth range r d = (10 , 10 
sorption should have little effect on radiation momentum 
transfer in dwarf galaxies and only becomes significant in 



ies (e.g. 


Lehnert & Heckman 


19961; Sanders & Mirabel 11996 


ICalzetti 


2001); however, see 


Krumholz & Thompson! (boi2 



for a counter-argument. In this maximal case, the optical 
depth Td is of order unity and could boost the momentum 
transfer by a factor of a ~few if the dust were to co uple with 
the gas fe.g. iMurrav. Quataert fe Thompsonll2005h . 

In these low metallicity or dust-depleted cases, Lyman a 
radiation can be trapped through resonant interactions, pro- 
ducing a possible signi ficant source of radiation pressure 
jHennev fc Arthurll998h . This process becomes insignificant 
when dust opacities reach unity, where dust reprocesses UV 
starlight into IR. Lyman a radiation accumulates in the res- 
onances, and eventually its pressure saturates, enhancing 
momentum transfer by a factor of 



/trap.Lya ~ 0.06 



(Td 



"d,MW 



He 



1/2 



(19) 



where R Bq is th e radius at where the gas pres sure equals radi- 
ation pressure l|Krumholz fc Matzneijl2009h . Thus we could 
be further underestimating radiation pressure as Lyman a 
trapping should become important below ~ O.O5Z , given 
similar dust properties as the MW. 



6 CONCLUSIONS 

We present results from three cosmological radiation hy- 
drodynamics simulations that demonstrate the role of ra- 
diation pressure in affecting star formation in dwarf galax- 
ies. We progressively add more physics to our simulations 
to understand the impact of each process. Our reference 
model considers primordial chemistry and Pop II and III 
star formation and feedback with a transition to a Salpeter 
IMF at 10 _4 Zq. In the next model, we add radiative cool- 
ing from metals. In our final and most realistic model, we 
add an lib-dissociating radiation background and momen- 
tum input from stellar radiation. Our treatment of momen- 
tum transfer from ionising radiation is accurately calculated 
from the optical depth given by our radiation transport mod- 
ule ENZO+MORAY down to our resolution limit of 1 comov- 
ing parsec. To our knowledge, this is the first cosmological 
galaxy simulation that has included the effects of radiation 
pressure that is computed from the radiative transfer equa- 
tion. 

The additional radiative cooling from metal species in- 
duced a runaway star formation event with a sSFR = 150 
Gyr -1 and SFR = 9.9 x 10~ 2 M yr~ x in our most massive 
galaxy with A/dm = 2.0 x 1O 8 M at z = 8.1. This strong 
burst creates enough ionising photons to reionise the entire 
1 Mpc 3 comoving domain; however, it creates a stellar popu- 
lation with a mean metallicity of 0.3 Z that does not follow 
the mass-metallicity relationship of dwarf galaxies. The lat- 
ter fact indicates that the galaxy is overcooling and overpro- 
ducing stars. Historically, this was avoided by using higher 
resolution and SNe feedback, but our simulations have both 
of these attributes. We find that radiation pressure plays a 
key role in regulating star formation in high-redshift dwarf 
galaxies. The highlights of our work on the effects of radia- 
tion pressure during the formation of a dwarf galaxy are as 
follows. 
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(i) Radiation pressure alone can reduce the SFR by a fac- 
tor of ~5 to 1.8 x 1CP 2 Mq yr" 1 , and a similar decrease is 
seen in the sSFR to 37 Gyr -1 . About 3 per cent of the ac- 
creted gas form stars in the buildup of the stellar population. 

(ii) The collision of dense shells that are driven by radia- 
tion pressure creates a turbulent medium with an rms veloc- 
ity equivalent to the circular velocity of the host DM halo. 
This added turbulence mixes the SNe ejecta thoroughly with 
the entire ISM. The turbulent pressure and reduced metal- 
licity prevents overcooling and the subsequent catastrophic 
collapse. 

(iii) Radiation pressure in dwarf galaxies is not sufficient 
to launch outflows, but it provides an additional impulse to 
drive dense gas away from massive stars and drive interstel- 
lar turbulence. SNe energy input is still the main mechanism 
to drive outflows in low-mass galaxies. 

(iv) Star formation regulation by radiation pressure 
creates a stellar population with a mean metallicity of 
10~ 21 Zq with a normal metallicity distribution with stan- 
dard deviation of 0.2 dex, agreeing with the mass-metallicity 
relationship of local dwarf galaxies. 

(v) We show that the majority of the momentum transfer 
occurs in dense shells at the edges of H n regions, i.e. it is 
not a local process. Furthermore, the timing of numerical 
radiative and SNe feedback processes is critical in capturing 
their effects. If the SNe feedback were to be released imme- 
diately into the ISM, densities will be artificially high and 
the SNe input will radiate too much of its energy. In nature, 
H II region dynamics will evacuate the immediate volume 
surrounding massive stars of dense and cold gas, and SN 
explosions occur in a warmer, more diffuse ISM. Our sim- 
ulations capture this temporal sequence, and the SN blast- 
waves do not overcool and are able to launch outflows from 
the galaxy. 

(vi) The column density of metals in 10 8 — 1O 9 M0 dwarf 
galaxies reach a maximum value of ~ 10 19 cm -2 . This may 
imply that that momentum transfer via dust grains is not 
important in low-mass early galaxies. 

We have shown that radiation pressure regulates star 
formation in dwarf galaxies in addition to photo-heating and 
SNe feedback. Its inclusion in galaxy formation simulations 
is key in forming realistic stellar populations and avoiding 
the overcooling problem also in high-resolution simulations 
that capture star forming regions with many computational 
elements. The impact of radiation pressure should continue 
to be important in larger galaxies, and we plan to apply 
our numerical methods to simulating such galaxies that are 
currently observed in the Hubble Ultra Deep Field. 
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